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The dimensional quantities referred to in this paper are 
given in U.S. Customary Units, in which the ntmierical 
calculations were made. The units are in brackats following 
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and parameters, the brackets are not included. 
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Chapter I 
INTRODUCTION 

On-board flight control and guidance is a subject which has 
had varying reception in different fields of Aerospace 
Engineering. In the area of unmanned missiles there has 
been extensive research, with many resulting applications, 
in developing on-board guidance systems, ais reported in the 
survey papers, Refs 1 and 2. These studies have encompassed 
many new optimal control and even differential gaming ideas 
(Ref 3 ) : in this field the on-board flight computer is an 
accepted and usually necessary part of the guidance system. 
While conventional homing and proportional navigation 
guidance laws are )le, and require minimal computation, 
more complex guidan.,^* schemes may be implemented on-board by 
the use of singular perturbation methodology, as in Ref 4. 

The willingness to apply state-of-the-art theoretical 
developments to manned aircraft is not as evident. This may 
be the result of a more conservative approach in applying 
new technology to machines which are responsible for 
peoples' lives, machines which are also extremely expensive, 
generally larger and more complex than many missiles. 
However one of the greatest obstacles may be the threatened 
removal of authority from the pilot: despite the existence 
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of sophisticated autopilots on many expensive aircraft, 
there is an aversion to total automation, particularly on 
the part of the pilot. As a result there is a significant 
gap between the flight-path optimization and differential 
gaming results which have been achieved in the last twenty 
years, and the their applications in on-board use. A part of 
this is due to the limited computational resources 
available, particularly on fighter and small general 
aviation aircraft, where weight and space are at a premium. 
Some of the latest developments relating to the latter case 
are given in Ref 5. On the other hand in the area of large 
transport aircraft the cost, weight and complexity of a 
small main frame computer is justified, but this has yet to 
be implemented. In civil aviation much research has been 
done in the area of trajectory optimization, with particular 
emphasis on efficient fuel usage and minimizing the direct 
operating cost. Attention has focused on the calculation of 
sub-optimal flight paths, using order- reduction to simplify 
the problem, as in Refs 6-9. Burrows (Ref. 6) used singular 
perturbations and order reduction to derive sub-optimal 
short and long haul trajectories, with on-board corrections 
to speed and energy errors based on expanding the 
performance index to second order, which he found to be more 
effective than simple linear feedback. Sorenson and Waters 
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(Ref. 7) used an assumed constant energy cruise (as did 
Erzberger and Lee, Ref. 8), and pointed oat that the on- 
board flight control needs to be coordinated with the ATC 
system, so that fuel saved during the flight is not wasted 
due to traffic congestion at the terminal area. Chakravarty 
and Vagners (Ref. 9) attempted -o provide justification for 
their state variable-selection through the use of non- 
diraensionalization. Transitions onto fuel-optimal climbs and 
descents are studied in Ref. 10, where they are used to 
derive a near-optimal feedback control law. Sub-optimal 
terminal guidance is examined by Erzberger, Ref. 11, for a 
fixed-wing aircraft, and by Baser, Ref. 12, for a tilt-rotor 
Aircraft. Optimal shipboard terminal guidance is studied in 
Refs 13-15. Despite the active interest and work, as 
described above, in this area the applications have lagged 
behind. A description, for example, is given in Ref. 16 of 
the DC-9-80 Digital Flight Guidance System; here the 
emphasis is on establishing reliability and safety criteria 
for the engine and flight control systems. It seems safe to 
say that in this area applications efforts have focused on 
feasibility and reliability rather than optimality. As 
mentioned earlier, the computational resources on a fighter 
aircraft are even more limited than on a transport, for 
obvious reasons of space and weight constraints. In contract 


with large transports there is a much greater range of 
applications for on-board optimal oontrel for fighV’^r 
aircraft. This io because a fighter can and often has to 
perform a much wider range of maneuvers (in terms of flight 
path angles and back angles for instance) as studied In Refs 
17-19. In many missions there is less, if any a priori 
knowledge of the flight path. Also it is often desirable for 
security to minimize the communication with the ground/ 
which eliminates the possibility of solving flight-control 
problems on the ground and relaying commands to the air. 

Witii this background it is the objective of this study to 
investigate on-board real-time flight control, with the 
intention of developing algorithms which are simple enough 
to be used in practice, for a variety of missions involving 
three-dimensional (3-D) flight. Ini-cially an approach is 
developed which is restricted to the intercept mission in 
symmerric flight, based on Ref. 20. Extensive computation is 
required on the ground prior to the mission but the ensuing 
on-board exploitation is extremely simple. The scheme takes 
advantage of the boundary-layer structure common in singular 
perturbations, studied in Ref. 21, arising with the 
multiple time scales appropriate to aircraft dynamics. 
Energy modelling of aircraft, as first examined in Refs. 
22-24 and extensively developed in Refs 25-27 is used as the 
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starting point for the analysis, 
nominal path is generated which 
cruise state. Feedback coefficients are found as functions 
of the remaining energy- to-go (dash energy less current 
energy), along the nominal path. These serve to generate 
transitions towards the nominal path, closed loop and to 
counter disturbances. In thir situation the guidance method 
is similar to the neighbouring-optimal guidance methods of 
Refs 28-32; these have been applied space shuttle re-entry 
problems, Refs. 33-35, and orbital transfer guidance. Refs. 
36-37. However there are two significant differences 
between this study and these references. In the present work 
the gain indexing is done in terms of the current energy; 
this avoids the problems encountered in estimating the index 
time, ai5 in the time- to-go or min-di stance methods. Also, 
for the extension to 3-d flight, families of reference paths 
are used instead of a single trajectory , with heading- to-go 
a.s the additional running variable. 


Urt-C.-v'- .j,, 

OF POCP. i 

In the symmetric case, a 
fairs into tlae dash or 
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1.1 PROBLEM FORMULATION 

The overall problem is to develop an on-board, real-time 
flight control system, which is near-optimal, for an 
aircraft flying an intercept mission, with arbitrary initial 
conditions. The equations of motion for a point-mass model 
of an aircraft can be written: 

E = .V(nT - D)/W 

h = VsinJT 

• 

y = (Lcos0 - Wcoslf)/mV 

X = Lsin^/mVcosy 

« 

X = VcosYcosx 
y = VcosYsinx 
m = -tjQ 

These equations embody the assumptions of thrust along the 
path, zero side-force, and flight over a flat earth with 
constant gravity. Also winds aloft are assumed to be zero, 
and the atmosheric properties standard. 
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1.2 SYiyiMSTRIC FLIGHT 

The first approach was to restrict the problem and simplify 
the model considerably, to reduce the analytical and 
computational burden, during the initial research and 
development of the guidance scheme. The restrictions in the 
problem, are the following: to consider only syimnetric 

flight, with fuel open, i.e. fuel optimization is not 
examined, which leads to maximum thrust in most maneuvers of 
practical interest. The target is assumed to be at a 
sufficient distance from the interceptor that a climb-dash 
is required: in other words a range-optimal climb to the 

dash point on the level flight envelope, blending into a 
steady- state dash. This sequence ends with a terminal 
transient, which is considered briefly in the next chapter. 
The time spent during the climb is assumed to be much 
smaller than the time spent at the dash state. The 
restriction in the aircraft model is that the variation in 
mass due to the fuel expenditure is ignored. Under these 
limitations, the ecfuations of motion are reduced to; 

• 

E = V(tiT - D)/W • 
h = VsinY 
y = (L - Wcosy)/mV 


X = Vcosy 


L 
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1.3 AERODYNAMIC MODELLING 

The aircraft which is used as an example to perform 
numerical calculations is a high-performance interceptor. 
The drag is modelled as a parabolic function of the control: 

^d = ^do ^ ^ ^'l^ 

The coefficients and K are functions of Mach Number: 

=do = “do*”’ 

and 

K = K(M) 

The thrust is a function of Mach Number and altitude; 

T = T(M,h) 

1 ! The way in which these three functions are represented is 

important in the computational work undertaken in this 
study. The reasons for this are discussed, . and the different 
methods which were used are described in Chapter 5 and 
Chapter 6. 


>• 


% if 
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Chapter II 

OPTIMAL CONTROL: REDUCED-ORDER MODELLING 

Reduced order modelling, based on time-scale separations 
observed in vehicle dynamics, is particularly attractive to 
the analyst in solving problems for lifting atmospheric 
flight. Numerical computations are simplified by the 

reduction in the system order and as a result the number of 
initial conditions which may have to be guessed or iterated 
upon is also reduced. Further, an improvement in the 

conditioning of the differential equations results from the 
confinement of the more unstable dynamics to boundary- layer 
^ corrections, which are relatively short in time. It has 

been appreciated since Kaiser's early work (Ref. 22) that the 
h and variables can be changed much more rapidly than the 
specific energy, E, which explains the introduction of this 
new variable. Also the energy can be thought of as a 'fast' 
variable in comparison to the range, at least in cases where 
the climb is a transient which fairs into a steady-state 
cruise or dash condition, i.e. when the time spent in the 
steady state is much greater than that spent on the climb, 
as assumed here. This leads to the reformulating the 
equations of motion, following the development of Ref. 25, 
with the inclusion of the interpolation parameters, on 
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the left hand sides of the differential equations for h and 
y, and E* on the left hand side of the differential equation 
for E; 

h = V SIN Z 
Z = (L - W COS 3T) / H V 
e^E = V(nT-D)/W 
X •= V COS Z 

To solve the problems of time-optimal ' control the 

variational Hamiltonian is formed: 

H = X_E + X, h + xJ + X X 
E h 4 X 

and the Maximum Principle (Refs. 38 and 39) is applied. The 

resulting Euler differential equations are: 

E- X. = -3H 
h 3h 

E^ Xy = -3H 
37 

« 

E^ Xp = 

3E 

t 

X„ = -3H 
3x 

The introduction of three separate time scales in the state 
system must conform to the requirement of the Tihonov theory 
(Ref. 40) that the ratio (e^/e^) 0 as e^ 0, as shown in 

Ref. 25. When both e^ and e' are equal to 1 the original 
point-mass model is recovered. 
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2.1 RECT I L I NEAR- MOT I ON MODEL 

The simplest model possible is obtained when both and c* 
are taken 0. By examination of the differential equations, 
the following consequences of these assumptions may be 
noted: 






_ 0 

h = 0 


y - 0 


II 

•Vo 


L = W 

= 0 

« 

E = 0 11 

T = D 


These equations embody the assumptions that the altitude, h, 
the path angle, Z, and the energy, E, can all be varied 
instantaneously in a control- like* fashion subject to the 
constraints. In this slow rectilinear-motion model the 
path- angle is, however, fixed at a value of zero, and the 
lift coefficient is chosen at ar.y energy/altitude 
combination so that the lift equals the weight. Further, 
the throttle is constrained so that the horizontal forces 
are balanced. The energy and altitude are chosen to minimize 
the Hamiltonian. This consists only of the range rate and 
the associated multiplier, which is constant because the 
Hamiltonian is not an explicit function of range in this, or 
any other modelling in this study. As a result the min-H 
operation leads to the high speed point on the level flight 
envelope. In the language *of singular perturbation theory 
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this is the zeroth-order 'outer solution', which the 
solutions from the other time scales must fair into 
assymptotically . The matching of different solutions and the 
composite generation are discused in a later section. The 
next time-scale is now examined. 


2.2 ENERGY- STATE MODELS 

The next level of order reduction is generally referred to 
in the literature as energy modelling. In this case is^ 
set to 1, and e^ to zero. Again the altitude and path- angle 
are assumed to be 'fast' and and 'control-like', but the 
energy change is analyzed and E assumes the role of a' 'slow' 
variable. Again the path- angle is fixed at zero, and the 
lift coefficient chosen so that the lift equals the weight; 
but the only remaining 'control-like' variable ( apart from 
the throttle, n) is the altitude; at any energy the altitude 
must be picked so as to minimize the Hamiltonian, which is 
now defined as; 

H = X^x + XgE = X^V Xj^E 

where the differential equation for E is given by: 

E = V(riT-D)/W 

Xp = -^ 

^ 3E 


The altitude which minimizes the Hamiltonian is therefore 
going to be determined, at any energy, by the relative 
values of X_ and X„ and their signs: their ratio determines 

the relative importance of range rate and energy rate, and 
their signs determine the sense of the optimization. For 
example, if X^ is small enough the altitude picked will 

correspond to the maximum possible instantaneous range rate 
possible at that energy, if X is negative. This is the 

lowest altitude (and highest speed) which is allowed" by the 
terrain limit, dynamic-pressure limit or Mach limit, 
whichever is greatest. On the other hand if the range 
multiplier is set to zero the altitude chosen will maximize 
the instantaneous excess power or energy rate, if X^, is 

negative. This special case is the so-called 'energy- 
climb', and is discussed in the following subsection. Note 
that if either multiplier is positive the rate of change of 
the associated state will tend to be minimized. 


2.3 ENERGY CLIMBS 

Of the' possible energy-state results the energy-climb is the 
simplest to calculate; as the Hamiltonian only contains one 
term, only one differential equation needs to be integrated 
assuming that X„ remains negative. The initial value of the 
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multiplier does not in general have to be determined; so 
long as it is negative the same path will result. Indeed if 
time histories are not required ncne of the differential 
equations need to be integrated at all: the altitude-energy 
path may be found simply by maximizing the level-flight 
energy rate at any energy. The energy climb for the aircraft 
studied is shown on Fig.lc It is interesting to note that 
this schedule shows multiple jumps in altitude, arising from 
realistic variations in the thrust data. This is somewhat 
different to other examples which have been examined, for 
example the F-4, where the altitude discontinuities in the 
energy-climb are primarily due to the transonic drag- rise 
(Ref 41). 


2.4 ENERGY- RANGE CLIMBS 

When the range multiplier, X , is not assumed to be zero, 

A 

i.e. 'energy-range climbs' are examined, the analysis and 
resulting computations are slightly more complex than the 
' energy- climb' discussed above. First of ail the X^, equation 

must be integrated, as the relative magnitude of X_ to X^ at 

X 

any time or energy is important in choosing the altitude. 
Secondly, as a result of this, the initial ratio of X^. to 

X , r" , must be carefully picked: different values of r“ 

X 
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will result in different paths with different terminal 
states. As the value of r" is increased from zero the 
resulting trajectories move downward in the flight envelope, 
with the terminal energy moving from the maximum energy, 
E , towards the dash energy, E, . At a certain value of r” 

= R® a path results which fairs gracefully into the dash- 
point. This is the range-optimal 'energy-range climb' which 
is desired and is shown in comparison to the energy climb 
found earlier in Fig. 2, with the level flight envelope also 
shown. Determining the correct value of r“ is an initial- 
value problem, but limited to only one dimension, and the 
usual one-dimensional search techniques, (i.e. golden- 
section, cubic and parabolic fits) may be employed. For 
values of r“ which are greater than R" the resulting 
trajectories are range-optimal for terminal energies which 
are lower than E^, over different time spans. These paths 

are characterized by a climb which approaches the dash 
point, a dash, and finally a terminal transient which takes 
the energy down to the desired level. This transient begins 
with an instantaneous dive to the maximum range rate (speed) 
at E^, as allowed by the terrain, dynamic-pressure, or Mach 

limit, whichever is the most severe restriction at the 
current energy level. In the case studied, no Mach limit 
and dynamic-pressure limits were applied; rather the thrust 
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data was faired off to limit the level-flight envelope from 
exceeding such limits, as explained in Chapter 6. As a 
result the terminal manoever takes the aircraft down to the 
terrain limit, (outside the flight envelope), where it 
remains, losing energy. This situation is unchanged until 
the energy is reached corresponding to the dash speed at the 
terrain limit. At this point the engine is switched off (Xg 

changes sign) and were speed brakes included in the model 
they would be applied: the instantaneous energy rate is made 
most negative. This sequence is shown in Fig. 3 for the 
aircraft being studied. For the case were Mach and dynamic- 
pressure limits are applied the equivalent maneuver is shown 
^ in Fig. 4. 

This process needs -some explanation: when the is less 

than E^, the aircraft must perform some terminal transient 

which loses energy in the most range-optlKnal way. There are 
two choices, or ways in which it can lose energy; at speeds 
below or speeds above the dash speed. Obviously the range- 
optimal strategy is to spend as much time in the latter 
region and as little in the former as is possible. This is 
done by switching off the engine when the speed drops below 
the dash speed, and if possible extending the drag brakes. 
The problem of the terminal-maneuver transient is not 
pursued here; it is of research interest. 
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2.5 METHOD OF MATCHED ASYMPTOTIC EXPANSIONS 
By the use of singular-perturbation theory, boundary- layer 
type corrections can be used to overcome the energy- 
modelling weaknesses, i.e. initial and final jumps in 
altitude, as in Refs. 25 and 42, and transonic or internal 
jumps, as in Ref. 41. While the altitude discontinuities 
are eliminated by expansion to the zeroth order, realistic 
path-angle values are obtained, in the Ref. 25 approach, 
only by continuing the expansion to the first order or 
higher. This is a nontrivial problem in the case where the 
altitude transitions occur at the beginning or the end of a 
trajectory, and is even more complex in the case of the 
internal jump. As a result, even the corrected energy model 
loses its attraction when realistic path-angles are required 
for onboard use as commands. A scheme for providing more 
realistic path-angle results in the zeroth order is explored 
in Ref 43 . 


2.6 CONCLUSIONS 

To conclude this chapter, some of the results of the 
reduced-order modelling are summarized below. 

First of all energy-state modelling, while attractive in its 
simplicity, is inappropriate and unsuitable for on-board 
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guidance use on its own, i.e. uncorrected, for the intercept 
mission contemplated. This is because it generates 
significant initial and terminal discontinuities in altitude 
and path-angle, which the aircraft is supposed to follow 
instantaneously. Secondly, multiple instantaneous jumps are 
also possible along along the optimal path, and lastly the 
path-angle is obtained as sero, in the usual approximation, 
which is a again a big disadvantage as the actual path- 
angles can be quite large. 

Corrections to the energy-state model which overcome 
these weaknesses are possible and have been demonstrated in 
the literature (Refs 41,42). Hci^ever this additional 
complexity is extremely unwelcome for on-board calculations 
due to limited storage and, more importantly, execution time 
available on-board; indeed solutions are not guaranteed due 
to the instabilities of the state-Euler system which need to 
be suppressed. In this context it is questionable whether 
this approach is in fact easier or quicker than solving the 
optimal control problem for the full system. 

However, certain ideas from the energy-state model are 
undeniably attractive. The solutions suggest a hierarchical 
structure of states in optimal control solutions. This is 
exhibited in the following way.* altitude and path-angle 
'command' values are determined by the current energy, and 
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in this sense the energy is the dominant state. If the 
current values h and JT do not coincide with these 
predetermined values a rapid transition can be made which 
brings them to their 'correct' values. These ideas form the 
basis of the guidance scheme which is presented in the next 
chapter. 
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Chapter III 
ON-BOARD GUIDANCE 

An alternative to using order reduction, suggested in Ref. 
20, which is simple enough to lend itself to onboard 
implemention is now developed, for the case of symmetric 
flight. The scheme has roots in the hiearchical structure 
of optimal-control solutions of the energy model, in which 
the specific energy is a relatively 'slow' variable and its 
value determines the control-like 'fast' variables, h and Z. 

3.1 NOMINAL PATH 

The phenomenon described above suggests that trajectories of 
the point-mass model funnel rapidly, (rather than 
instantaneously as in the energy model), into the vicinity 
of a single path, which leads to the dash-point. The idea 
pursued in this thesis, and Ref 44, based upon Ref 20., is 
to determine this 'skeletal path' for the point-mass model, 
for as wide a range of energies as possible. This is the 
nominal, or reference trajectory and the altitude and path- 
angle histories are recorded as functions of the energy or 
energy-to-go, rather than time or time-to-go, as is common 
in other neighboring optimal guidance schemes (Refs 28-37). 
The advantage of this approach in an on-board context is 
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that approximations to the final time are not necessary, and 
implementation of the scheme is greatly simplified as a 
result. 

3.2 FE EDBACK LAW 

The next step is to generate a neighboring-optimal feedback 
guidance law which will control the aircraft so as to follow 
a neighbor of the nominal optimal path. There are two basic 
reasons for doing this. First of all the reference path is 
of little use open loop: even if the aircraft is at any time 
on the reference path, the control commands which are stored 
along this trajectory will be insufficient to keep the 
aircraft close to it. This is because disturbances and 
errors inevitably arise both in the actual flight (i.e. 
variable winds etc) and in representing the control history 
using a cubic spline (Ref 45) Secondly, even if this first 
problem could be ignored, the reference path is of little, 
if any, use when the aircraft has initial conditions which 
are far removed from the nominal: for instance if the 
aircraft is initially loitering at high altitudes and 
subsonic speeds, on combat patrol, for example. Linear- 
feedback coefficients are proposed to generate the necessary 
transients to bring the aircraft to the neighborhood of the 
nominal optimal and stabilize the subsequent path. The 
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guidance law is a linear feedback control based on the 
difference between the nominal and actual altitude and path- 
angle values. 

3.3 FEEDBACK COEFFICIENTS 

The feedback coefficients, which correspond to minimizing 
the second-variational approximation to the performance 
index, as in Refs. 28-37, are found by perturbing the 
altitude and path-angle separately from their nominal values 
along the reference trajectory. The optimal-control problem 
is re-solved and the partial derivative of the control with 
respect to the states {at fixed energy) is estimated by 
difference quotient approximation. The partial derivatives 
which are mentioned here are the variations in the 
parameters of an initial value problem; they should not be 
confused with the variations of the control along the 
trajectory. They are defined for an arbitrary value of 
energy = E^ in the following way: 

let C^*{t) be the control which takes the aircraft from an 

initial point at low energy, E®, (altitude and path-angle 
zero), along the nominal path up to the dash point on the 
level flight envelope, while optimizing range; the resultant 
state time histories are given by 


h*(t), y*(t), and E*(t) 
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let the energy of the aircraft reach the value while 

travelling along the nominal path, at a time t^: 

E‘ = E*(t*) 

Then at E^ the 'correct' altitude and path-angle are given 


by h*{t^) and 

y*(tM- To 

find 

the 

altitude 

feedback 

coefficient at 

this energy 

level 

the 

procedure 

is as 

follows : 






find the range- 

“Optimal path 

which 

has 

the same 

terminal 


conditions, and terminal time as before but use the nominal 
state at t^ -as- the initial conditions, with a perturbation. 
Ah, introduced in the initial altitude; 

r(0) = r*(tM 

E(0) = 

h(0) = h* (t‘ ) + Ah 

The solution of this problem results in a new control time 
history, altitude feedback coefficient is 

found by the following secant approximation; 
ah Ah 
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3.4 ON-BOARD USE 

The commands to the autopilot are tak^n from the nominal 

path with linear corrections for the variation of the 
altitude and path-angle from their nominal values. On-board 
use requires only the storage of the states (h and T) , 
control (lift coefficient or load factor), and the two 
feedback coefficients, each as functions . of energy, or 
energy- to-go. The feedback guidance law with the appropriate 
functional dependencies are shown below: 

C, = C,*(E) + |C, (E) ( h - h*(E)) + 3C, (E) ( r - y*(E)) 

' ' 3h' ay' 

To summarize the only variables required to be stored on- 
board in symmetric problem are : 

C^*(E) 

h*(E) 

y*(E) 

3C,(E) 

ah' 


(E) 

ay ‘ 


Chapter IV 

OPTIMAL SOLUTIONS FOR THE POINT-MASS MODEL 

A requirement of the proposed idea is a large number of 
optimal-control solutions to the point-mass-modelled 
problem. Optimal control solutions can be found in many 
different ways. They can be found by the use of direct 
methods, such as gradient methods, where the control history 
is parameterized in sectionally-linear or spline 
approximation and the terminal conditions are met by either 
penalty or projection techniques. Alternatively, the problem 
can be resolved into a two-point boundary value problem, 
with split boundary conditions. Half are known at the 
initial time and the other half at the final time. This can 
be solved by the use of indirect methos such as simple or 
multiple shooting (Refs. 22,23). To solve the problem of 
time-optimal control the variational Hamiltonian is formed: 

H = * Xjjh + 

and the Maximum Principle (Refs 38,39) is applied. 
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The rsulting Euler differential equations are : 

L = -8H 
^ iE 

Xh = -3J3 
3h 

iy = -8H 

ay 

x^ = -9H 

ax 

The lift and the throttle setting must be chosen to minimize 
the Hamiltonian, which requires that: 

m = 0 
aC"! 

and 

(X^ < 0) 


4.1 


METHOD OF SOLUTION 


Euler solutions were found in the present work by the method 
of multiple shooting, using the algorithm and computer 
program of Refs. 33,48 kindly supplied by DFVLR, 
Oberpfaf fenhofen. West Germany. In this method, the 
interval of integration is broken up into many subintervais. 
This is preferable to 'simple shooting', where the initial- 
value problem is attempted directly, as optimization 
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attempted directly 


problems of lifting atmospheric flight are ill-conditioned, 
the state-Euler system being violently unstable. 
Partitioning the time interval has the effect of suppressing 
error growth. This method was used primarily for reasons of 
accuracy. This need arises, for example, in the calculation 
of the. feedback gains, found by the difference of the 
control at the beginning of two optimal solutions. 
Typically to find the gains to 5 fii^res the control must be 
known to about 8 figures. The* multiple-shooting method has 
greater' accuracy than the other methods available, and 
although it is often difficult to generate the initial 
reference trajectory, the subsequent calculation of the 
feedback gains is relatively easy as the method has good 
convergence properties in the vicinity of a solution. 
Further discussion on these topics is found in Chapter 8. 





Chapter V 

INITIAL EXPOSURE TO OPTSOL 

The first use of the multiple shooting program OPTSOL 
obtained from DFVLR was to solve a very simple optimal 
control problem. This, taken from Bryson and Ho (Ref. 49) 
page 121, is similar to the brachistochrone, and was solved 
numerically both with and without a constrained arc, to test 
the user- supplied software required for the program. 
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5.1 AIRCRAFT DATA MANIPULATION 


The program 

OPTSOL 

had been 

brought 

to 

VPI&SU with 

subroutines 

already 

created to 

enable 

the 

solution of 


aircraft flight- mechanics problems and, rather than try to 
start from the beginning, attempts were made to use the 
existing computational tools, at least until familiarity had 
been gained with the program. In particular, the data which 
was used to model the aircraft under study was extensively 
modified so that the integration subroutine in OPTSOL, known 
as DIFSYS, was able to function. This proved to be a 
problem, as DIFSYS, as received, was extremely sensitive to 
the degree of smoothness of the right hand sides of the 
differential equations. In fact if discontinuities are 
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encountered in any derivative up to the eighth, the stepsise 
of integration shrinks to zero. As all data of the point- 
mass model had been represented by cubic splines and spline 
lattices to facilitate interpolation, considerable effort 
was spent on the generation of an analytical representation 
which would reproduce both the values and the shapes of the 
data with consistency. This had been done at DFVLR by using 
polynomial expressions, and this method was examined for the 
aircraft data on hand and abandoned. While a polynomial of 
sufficiently high order will fit any number of consistent 
data points exactly, there is an increasing distortion of 
shape with increasing order of polynomial. In fact even low 
order polynomils did not match the data ' at all well. The 
approach taken was to use a combination of polynomials, 
exponentials and arctangent functions to accomplish this. In 
the case of the single valued functions, i.e. K(M), 

this was not too difficult. The arctangent functions can be 
used as' soft’ switches, separating different portions of the 
data, which can be represented by a simple function locally 
(i.e. by a straight line or a parabola). However in the case 
of multivariable functions such as thrust and fuel flow this 
is definitely a nontrivial problem (however only thrust was 
attempted) . In the case of thrust, the representation was 
achieved by fitting against Mach number, using coefficients 


which were functions of altitude. 19 variables were 
optimized using a conjugate gradient process which minimized 
the sum of the square of the errors at the grid points . The 
functions developed for Thrust, shown in 

Appendix A, and the aerodynamic data are shown graphically 
in Figs. 5 and 6. 

After construction of the smooth data, the flight 
envelope was calculated and drawn (Fig 7). As in the case of 
some high-performance jet-fighter aircraft the envelope 
turns out not to be performance limited, i.e. the level 
flight maximum sustainable speed is much higher than the 
Mach limit. In this case M=2.4 is the Mach limit and the 
high speed point occurred at roughly M=3.0. It should be 
mentioned that aerodynamic and thrust data are not actually 
available for M= 2.4 and the flight envelope found by 
extrapolation is essentially a conjecture. The important 
thing is that the excess power at level flight is greater 
than zero for a range of altitudes along the Mach limit, for 
which both thrust and aerodynamic data are reliable. This 
problem, which in general requires treatment of state- 
inequality constraints, was dealt with in the following way: 
the Thrust was faired off sharply against Mach Number, near 
the Mach limit so that the flight envelope no longer 
exceeded it. This was done, by multiplying the thrust by a 
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switching arctangent function which rapidly (but smoothly) 
brought the thrust to zero while leaving it unaffected 
elsewhere. The dynamic -pres sure limit was treated in the 
same way. The analytical formulation for these two limits 
are included in Appendix A. The flight envelope with the 
Mach-number limit is shown in Fig. 8; the effect of both of 
the limits is shown in Fig. 9. 

5.2 INITIAL FLIGHT-MECHANICS PROBLEM 

Once the dataset had been finalized, OPTSOL was used to 
generate some optimal trajectories for a simple atmospheric 
flight problem: maximize final speed, from a given initial 
state, with final path angle zero and final altitude free. 
This was was solved for several different time intervals, 
using simple shooting (initially), and also multiple 
shooting, to gain familiarity with the use of multiple 
shooting and to investigate the methods of finding families 
of trajectories , for instance by time stretching. The time- 
histories for a family of four different trajectories are 
shown in Fig. 10-12. These are, respectively, speed, path- 
angle, and altitude plots. 
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3.3 FIRST TRAJECTORIES TO TOT DASH POINT 

The next step was to attempt to find paths which went to the 
high speed point, over a fixed time interval and to try to 
decrease the initial energy while lengthening the overall 
flight time. This was done by starting at an altitude and 
speed combination , (path-angle zero), just below the dash 
point, guessing the values of the costates. A total 
integration time of 5 seconds was used, and as can be 
imagined, the first guess was far from the targeted fin»i 
conditions; however by requiring OPTSOL to satisfy boundary 
conditions by successive proximity rather than in one jump, 
a trajectory which reached the specified altitude and path 
angle combination was found. However, it was not possible to 
to get the final speed to the desired value in the 5 second 
interval, because the time was not long enough to reach it. 
To achieve the desired final speed and to observe the manner 
in which the system approaches the equilibrium point (the 
possibility of an oscillatory solution near the high speed 
point , analogous to oscillatory cruise solutions was 
considered a possibility), attempts were made to lengthen 
the time of integration, by stretching the sub- intervals in 
the multiple- shooting scheme. Initially it was found to be 
very difficult to e.xtend the trajectory at all - OPTSOL 
would not converge for even extremely small increases in the 


final time. Eventually the interval was increased to 6 
seconds. The final speed also increased but still did not 
reach the value at the dash point. It became virtually 
impossible to increase the final time any further due to 
numerical integration difficulties. For this reason and 
computational expense, the approach was reassessed at this 
point. 


5.4 ' eigenvalue ANALYSIS 

The system was lineariaed about the high speed point to 
examine the dynamics of the system in the vicinity of the 
equilibrium point . The analysis revealed that the 
stablility eigenvalues were all placed along the real axis. 
At first the absence of complex roots akin to phugoid 
oscillation suggested that the linearization had been 
incorrect. After this had been checked and rechecked, the 
analysis was repeated at a point removed from the vicinity 
of the sharp arctangent functions which had been used to 
limit the flight envelope, it was conjectured that the 
switching functions may have introduced large gradients 
affecting the dynamics of the closed- loop system . The 
throttle coefficient was reduced to 0.68, reducing the speed 
of the dash point by about 100 ft/sec, well away from the 
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arctangent switch region, and the linearized analysis was 
repeated. The eigenvalues were found to have both real and 
imaginary parts, as expected, showing that the steps taken 
to limit the flight envelope had engendered significant 
effect on the dynamics of the state-Euler sytem. The s- 
plane positions of the two cases are shown in Fig. 13 and 

14 - 

5.5 BACKWARDS INTEGRATION OF STABLE EIGENVECTORS 

It was thought that a useful starting trajectory could be 
found by using the stable eigenvectors of the linearized 

' sytem. If the equilibrium state is disturbed in proportion 

to a stable eigenvector the disturbance will die out in the 
linear case and should fair in towards the equilibrium 
point, for some finite time at least, in the nonlinear case, 
if the disturbance is small enough. So if such a trajectory 
is integrated backwards in time (using the full nonlinear 
system) a series of points will be generated which will fair 
in towards the dash point, at least for some time. Only one 
of the three eigenvectors approached the dash point from the 
desired direction, i.e. from points lower in altitude and 
slower in speed. This was integrated for 22 seconds and used 
as an initial guess for -OFTSOL. The path-angle at the 
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initial time was non-zero and attempts were made to reduce 
it to zero. Again convergence troubles were encountered: 
OPTSOL could not tolerate large changes in the initial 
values and the effort was finally abandoned. Apart from the 
cost of computing and poor convergence behaviour, the system 
also displayed an alarming instablity to small changes: on 
occasions the speed in the final seconds dropped from its 
maximum value (about 2300 ft/sec) to 1 ft/sec. 


5.6 CONCLUSIONS 

It was concluded that the thrust-tailoring approach taken to 
make the problem easier had instead probably made it worse. 
The integration subroutine DIFSYS is very sensitive to small 
changes in derivatives of the right hand sides. By using a 
multiplicity of sharp arctangent functions the computational 
burden became large, as every time DIFSYS encountered an 
arctangent transition the stepsize of integration 
automatically became very small , increasing the computer 
time required. Further it was evident the system was overly 
sensitive to small changes in initial values. As a result it 
was decided to use a simpler integration subroutine and to 
return to splined data. 




Chapter VI 

MODIFICATIONS TO OPTSOL 

The first step to modify the operation of the program OPTSOL 
was to change the integration routine. The variable step, 
eighth order Runge-Kutta package DIFSYS seemed to be a 
primary scource of the numerical difficulties and 
computational expense experienced in the early use of 
OPTSOL. It was removed in favor of a much simpler fixed 
step-size fourth order Runge-Kutta-Gill subroutine. 


6.1 SPLINED AIRCRAFT DATA 

This substitution enabled the use of cubic splines and 
spline lattices of Ref (45) for representation. of the 
aircraft thrust and aerodynamic data. The problem of the 
Mach-limit violation was handled by fairing off the thrust 
data gently over four tenths of a Mach Number and increasing 
the drag by adding more missiles. The aerodynamic and 
thrust data are included in tables 1-4. The new flight 
envelope was calculated and is shown in Fig. IS. The 
coordinates of the dash point were found by a Newton 
iteration applied to the usual necessary conditions. 
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6.2 FAMILY OF TRAJECTORIES TO THE DASH POINT 


The new data were used to calculate an ' ene rgy- climb' 
schedule(Ref .25) ; this was used as a guide for guesses of 
initial altitude, energy and trajectory time combinations. A 
thirty-panel division of the trajectory was employed to find 
trajectories starting at lower altitudes, over longer times. 
This procedure was successful in finding optimal-range 
histories starting from an initial energy of 30,000 ft. 
After this point it became difficult and expensive to 
progress any further down in altitude and energy. It was 
thought that a smaller stepsize might be necessary to 
evaluate partials with sufficient accuracy for the method to 
converge. However this did not improve matters 
significantly. But when the program was brought to Langley 
Research Center the situation improved. The CDC computer has 
a word-length which is approximately double that of the IBM 
370, so with double precision at Langley about 28 decimal 
digits were obtainable compared to 14 or 15 digits at VPI. 
This had a significant effect on the program's operation. 
Much smaller stepsizes were used to evaluate the Jacobian 
without a penalty in round-off error, and it is conjectured 
that the resulting improvement in the accuracy of the 
Jacobian helped the convergence of OPTSOL. The trajectory 
extension continued until zero altitude was reached over a 
trajectory of 282 seconds. 


Chapter VII 

OPTIMAL- REFERENCE- PATH C ALCUL AT IONS 

The first objective is to generate a reference optimal path 
using point-mass-model dynamics, over the widest possible 
energy range. In the climb-dash problem, the highest energy 
of interest corresponds to that of the high-speed point on 
the aircraft envelope, the dash ’outer' solution. The 
lowest energy corresponds to the trajectory which just 
kisses the terrain limit, i.e. below this energy, optimal 
solutions which start at zero altitude would dive below the 
terrain limit if it were absent. This lower energy is found 
by examining the initial load factor of a family starting 
from level flight at the terrain limit altitude: when the 

initial load factor is unity the lower energy is determined. 
This is shown in Fig. 16, where the initial load factor is 
plotted for several different initial energies. 


7.1 FINAL LOAD FACTOR 

Once the energy had been found which pulled off the ground 
with an initial load factor of 1, the effect of the flight 
time was investigated. To satisfy the final conditions in a 
finite time rec[uires that the aircraft perform some 
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maneuvering near the terminal energy: the longer the time 

allowed to approach the equilibrium point, the more gradual 

the approach should be. The effecr of flight time on the 

final load factor was studied (for the same inital and final 

conditions) and results are shown in Fig. 17 . This clearly 

demonstrates how the optimal path tends to fair in 

asymptotically ?.s the flight time is increased. The load 

factor dropped to 1.001 after the flight time had been 

increased to 360 seconds. This time was chosen for the 

nominal path adopted in guidance- scheme development, and the 

altitude and path-angle (state variables) as well as the 

lift-coefficient (control variable) have been splined as a 

function of the energy. The load factor is shown in Fig. 18 

, drawn against energy, showing the grid points used in the 

spline. Fig. 19 - 22 show the energy histories for path- 

angle, altitude, load factor and lift coefficient 

respectively for t = 360 secs. The other paths from the 

max 

same initial energy, but over longer times, showed identical 
state and control energy histories over almost all the 
energy range. However at the terminal energies the effect of 
different flight times is most evident. Comparisons of the 
trajectories which result for different flight times are 
shown in Fig. 23 - 26 for path-angle, altitude, load factor 
and lift coefficient respectively. These variables are 
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plotted versus energy for the last 2000 ft of energy, for 

„ = 300 and t = 360 seconds. The dramatic effect that 
max max 

the flight time has on the final state and control behaviour 
is obvious from these pictures. 


7.2 ONE PANEL INTEGRATION 

After each converged solution was obtained a trajectory was 
performed for the entire time, from the initial conditions. 
At higher energies and over shorter times this would 
ordinarily generate final states which were close to those 
specified in OPTSOL, but owing to the error propagation of 
the mismatched paths at each grid point, there is a 
difference between a one-panel integration and a 30-panel 
integration. However at energies with zero initial altitude 
the error propagation was such the final conditions were not 
nearly met. After about 150 to 200 seconds the 
instabilities in the state-Euler system would produce 
extreme results. This raised the question as to whether the 
solution generated by OPTSOL is optimal or even near 
optimal. To this end the number of panels was reduced first 
to 10, then to 6. Attempts to drive the number smaller than 
this were not succesful as it appeared that the computer was 
' rvinning out of digits', despite the fact that 28 were being 
used. However the difference between the solution for 6 


panels and for 30 panels lies beyond the 9th digit and so it 
was assumed that no benefit would be gained by trying to 
reduce the number of panels. 

7.3 ENERGY-MODEL POINT-MASS-MODEL COMPARISONS 
Having established the nominal optimal path which takes the 
aircraft up to the dash point, it is of interest to stop and 
consider the two different models which have been used to 
study the problem, in particular it is of interest to 
compare the two different paths which climb up to the high- 
speed point. These are shown in the h-v plane in Fig 27, 
surrounded by the level-flight envelope. The energy-range- 
climb model is indeed close to the point-mass model 
particularly at higher energies. 




Chapter VIII 

FEEDBACK COEFFICIENTS - CALCULATIONS 

This chapter describes the numerical work done to evaluate 
and represent the feedback coefficients used in the guidance 
law for the case of symmetric flight. In this case the 
coefficient*’ are the partials of the lift coefficient with 
respect to the altitude and path-angle, at a fixed energy. 


8.1 METHOD QE EVALUATION 

The calculation of the variation in the control due to 
errors in the altitude and path-angle is treated as an 
initial-value problem, and has been extensively discussed in 
Chapters To improve the accuracy of the feedback 
coefficients, each one was evaluated twice, by introducing 
positive and negative perturbations, and taking the average 
of the two difference-quotient values. This method also 
allowed the determination of the optimal size of disturbance 
(in terms of the resulting accuracy) by varying the size of 
the disturbance , examining the degree of agreement between 
the two values untill the 'best' stepsize has been found for 
both altitude and path- angle. While it is true that the 
optimal stepsize will in general vary along the reference 
path, it was found that this change was negligible and one 
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value was effective in evaluating the entire range for 
either coefficient. As the stepsize is reduced the errrors 
due to nonlinearities shrink, but those due to a finite 
word-length grow: hence a compromise defines the optimal 

disturbance. It ha 3 been noted that a multiple shooting 
method such as OPTSOL is well suited to these kinds of 
calculations: although it was an arduous task to establish 
the nominal path, once this had been achieved, the 
neighboring solutions were found rapidly (within 3 or 4 
iterations) and with highaccuracy. This last point is 
important, as the use of numerical differentiation of the 
initial control to find the feedback gains required high 
precision the control. Typically it was found that 8-9 
decimal digits of information were required for 4-6 figure 
accuracy in the gains. 

8.2 PILOT SCHEME 

Feedback coefficients were initially found over a small 
range of energies, to evaluate the usefulness of the scheme 
before committing the computational resources needed for the 
full-scale operation. The last fifth of the energy range was 
chosen for this purpose as the integration times are the 
shortest and this minimises the cpu time rquired to find 
optimal control solutions. The energies and corresponding 


times were taken from the reference trajectory (of 360 
seconds) in the following manner; the total energy change 
was divided into twenty. The reference path was then 
integrated again and whenever the energy at the end of an 
integration step exceeded an integer number of divisions of 
the total energy change, the time and energy were recorded. 
The times and energies for the pilot run are shown in 
table?. The disturbance sizes were varied so as to maximize 
the agreement in between the two values obtained for each 
coefficient. The optimal perturbation -in- altitude was found 
to be 0.05 feet; in path- angle it was found to be 0.0000001 
radians. Agreement between the values of both of the 
coefficients was found to vary in between 4 and 6 figures. 

In addition to the energy levels already chosen for feedback 
coefficient evaluation, it was necessary to find values 
close to the final energy as well. This is because spline 
representations are very unreliable when used to extrapolate 
data. The energy at the beginning of the last panel in the 
multiple-shooting method, i.e. at 348 seconds, was chosen as 
the upper limit for this purpose. The gains at this energy, 
which is just 0.11 feet below the maximum value, turn out to 
be an order of magnitude larger than the gains at lower 
energies. This sensitivity of neighboring-optimal-guidance 
schemes close to the terminal state has been noted in the 


literature (Ref s. (28-37) ) . It is worth commenting , 
however, that the apparent unboundedness in the gains near 
the final state could have been a result of the method by 
which they were calculated; it is quite possible that a 
finite integration time, which is shorter as the terminal 
state is approached, was responsible. In other words if a 
longer time of .integration had been allowed for the paths 
which were close the final state a different behaviour might 
have been observed. However this effect is highly local, and 
due to limitations of time and money this, topic was not 
pursued. Any actual implementation of the scheme would of 
course have to take this into account, possibly by setting 
an upper limit on the magnitudes of the ^^ains, to avoid 
control saturation with small errors. To examine the 
transition in the feedback coefficients near the terminal 
state, the analysis was repeated for 3 more energies close 
to the final time, at 336,324, and 300 seconds. This is an 
inexpensive set of calculations as the integration times are 
extremely short. Also the coefficients were evaluated at the 
energy corresponding to the trajectory time of 188.7 
seconds, as it was felt that they were needed for accurate 
spline representation. 

The next problem was to spline the coefficients as functions 
of the energy-to-go. Difficulties were encountered at first 



when the splitting was attempted. Cubic splines are not 
suited in general to represent functions where large 
variations in the gradient exist. In this case rhe gradient 
changes by six orders of magnitude in the vicinity of the 
end-point, resulting in large extraneous oscillations 
appearing throughout the spline representation, which render 
the interpolation useless. One way (not very satisfactory) 
is to ignore the spurious points which are causing the 
trouble. This was done in this case, and the plots of the 
coefficients are shown Figs 28,39. 

To overcome these difficulties the splines-under-tension of 
Ref. SO were used. These are similar in character to the 
cubic splines of Ref 45 which had been used so far; the 
additional feature of the splines-under-tension package is 
the ability to minamise spurious wiggles near regions of 
rapidly changing gradient by the use of a tension factor, o. 
By increasing o ■ the anomalies can be reduced but not 
eliminated, at least in the vicinity of the end point. The 
problem is that as the tension factor is increased the 
oscillations near the end point die down but the rest of the 
representation becomes essentially polygonal, i.e. linear 
interpolation between the data points. 


8.3 LOGARITHMIC SPLINING 


It became apparent that the normal or ordinary method of 
splining was inadequate and a different approach was needed 
to continue. Essentially this is a boundary- layer type 
problem; there is a region where the coefficients vary 
rapidly. It seemed to be appropriate to separate the two 
regions and, using different methods, spline each one 
separately. The only requirement would be that the two 
representations fair into each othtvr smoothly. One 
possibility is to use the normal splines in the 'outer' 
region, and spline the terminal coefficients in terms of the 
logarithm of the energy- to-go, matching the slopes at the 
junction between the two regions . (Another possibility is to 
use the inverse of energy-to-go in the terminal region, but 
this was not used for reasons as the large variations in the 
gradients, which are the roots of this problem, still 
exist) . The logarithmic method was used to spline the 
coefficients for the range of energies considered in this 
pilot section. The results are shown in Figs 30 and 31. 
These show the gains using 10 gridpoints for interpolation. 
These show a dramatic improvement over the previous attempts 
to spline the data: these earlier efforts had been so bad 
that they would only be visible on the same graphs as a 
series of vertical lines passing through the gridpoints. It 
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was considered likely that with a few additional points the 
small remaining anomalies would be eliminated. An 
additional 16 points were evaluated in the vicinity of these 
outstanding 'wiggles' and finally a usable representation 
was generated, shown in Figs 32 and 33, as functions of 
energy. They are shown as functions of the logarithm of 
energy-to-go in Figs 34,35. 

When the decision was made to carry on and evaluate the 
coefficients over the rest of the energy range, the same 
method was used to spline the data: the logarithm of the 
energy-to-go was used, and there was no need to go to a 
boundary layer type of approximation after all. The 
coefficients as they were represented over the entire energy 
range are shown as functions of the energy in Figs 36 and 
37. The corresponding plots versus the logarithm of energy- 
to-go are shown in Figs 38 and 39. 


Chapter IX 

SIMULATION AND TESTING 

Following the satisfactory splining of the nominal states, 
controls and feedback coefficients as functions of the 
energy-to-go, the guidance scheme was tested by running a 
simulation of the point-mass-model, using the feedback law, 
and comparing the resulting trajectory with an Euler 
solution which started from the same initial conditions. 
Before the entire range of feedback coefficients had been 
worked out a pilot scheme tested out the idea on a small 
range of energy near the dash-point. This test was performed 
with an initial disturbance of 1000 ft; the trajectory which 
resulted from the guidance law is compared with the Euler 
solution from the same initial conditions and the nominal 
path in Fig 40 where the altitude is plotted as a function 
of energy. The guidance law is so close to the optimal path 
from the same starting point that it is almost impossible to 
discern the difference between them on this Figure. The 
difference in altitude between the two is shown as a 
function of time in Fig 41 it can be seen that the 
difference is always less than 11 ft. With zero disturbance 
the autopilot was able to follow the nominal path more than 
satisfactorily, over the entire range of energies, despite 


the inevitable errors which arise in the spline 
representations. Tests were perfc rmed with the initial 
altitude disturbed from that of the nominal path at 
different energies by 1000, 5000, 10000, and 15000 feet 
above and by 5000, 10000 feet below the nominal path. The 
resulting trajectories are shown in Figs. 42-46. These show 
that the feedback law follows, the optimal solution closely, 
even when the initial disurbance is far outside of the range 
of linearity of the feedback gains. The cost was calculated 
for the situation with an initial altitude of 15000 ft above 
the nominal value, at the point where the two trajectories 


faired 

into 

the 

dash 

point. 

The 

difference 

between the 

ranges 

was 

less 

than 

600 ft, 

an 

extremely 

small number 

considering 

that 

the dash speed 

is 

2400 ft/sec. 




Chapter X 

EXTENSION TO 3-D FLIGHT 

This chapter descibes the work done to extend the analysis 
to three dimensional flight, and suggests what direction 
future efforts might take. 

10.1 CROSS RANGE CONSIDERATIONS 

The problem of extending the analysis to 3-D flight is now 
considered. The state system is augmented to include y, the 
cross range, and x, the heading angle. The addition of the 
corresponding multipliers to the full state-Euler system 
raises the order of the problem to twelve. For the 
intercept problem the final value of y must be zero; the 
value of the final heading, relative to the initial heading, 
must either be calculated on-board, or be supplied by the 
GCI. This will in general vary, for a maneuvering target, 
and the value stored on-board must be periodically or 
continuously updated. 

The boundary condition on y leads to a dependence of the 
optimal solution on the cross range: for the same heading- 
to-go and energy-to-go there will exist many different 
possible values of y. As a result, if this formulation is 
used, cross range-to-go is an additional running variable; 
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this increases the order of the nominal paths required, 
which means a large increase in the computations on the 
ground, as well as an increase in the storage requirements 
on-board. 

To get around this situation it is proposed to avoid 
using an additional running variable by letting the final 
value of y be free: this can accounted for in the 
computation of the final heading needed for intercept, as 
specified by the on-board flight computer or the GCI. The 
intercept paths which result from the two different methods 
are compared in Figs 47 and 48, for a target which is 
initially far away from the interceptor. 

10.2 COMPUTATIONAL CONSIDERATIONS 

The first approach considered to generate, a family of paths 
to the dash point was to use the symmetric flight reference 
path as a starting point for the augmented system, and 
introduce a small heading- to-go at the initial time. The 
argument for doing this is that for very small headings the 
state-Euler system should not be changed very much: the 
paths are close to each other. However this method is only 
useful for a small number of combinations of heading- to-go 
and energy-to-go . This is because the turning rate at the 
energy at which the aircraft lifts off the ground is so high 
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that all the heading- to-go disappears in a short time, and 
over a very small energy range. In general a method must be 
found which generates the part of the family of reference 
paths which combines moderate and large headings-to-go and 
moderate to small energies to go. The difficulty lies in 
knowing what initial conditions to pick for the altitude and 
the path-angle; when the aircraft Is lifting off the ground, 
these variables are specified, but in the general case, 
starting from an arbitrary energy-to-go and heading- to-go 
combination, the selection is a problem. Letting them be 
free is not acceptable as it can lead to an initial lift 
coefficient of zero (i.e. in the symmetric case): the 
optimization algorithm takes advantage of the freedom to 
choose the initial conditions in a way which maximizes the 
short term benefit. This does not fit in with the concept of 
a nominal refence path, where the altitude and path-angle 
are the same at the same combination of energy and heading- 
to-go . 

The solution that is reccomended is to use the altitude that 
comes out of the energy- turn model, as in Ref 25. Here the 
heading is assumed to be a 'slow' variable, and has the same 
status as energy. However, instead of having to choose one 
variable, (such as the ratio of the initial energy 


multiplier to the range multiplier, as in Chapter 2), the 


initial heading multiplier must also be iterated upon. This 
is done using a Davidon-Fletcher-Powell algorithm, to find 
the path which fairs into the dash point with zero heading. 
An example of such a path over a small range of energy and 
heading-to-go is shown in Fig 49, where the heading is shown 
against energy, and in Fig 50, where the heading vs time 
plot for the same initial conditions is shown. 


10.3 - SELECTION OF Tlffi INITIAL PATH- ANGLE 

The energy- State model produces altitude predictions which 
are fairly accurate as a function of the current energy, 
(away from altitude jumps), as can be seen Fig where the 
Euler solution to the climb-dash is compared to the energy- 
range solution. However the same can not be said for the 
path-angle, which is predicted to be zero along the path. As 
a result a modification is considered, (Ref 43), which 
produces realistic values along the path. The difference 
lies in the selection of the fast and slow variables: if 
altitude is chosen, zero path-angle results, if velocity is 
chosen, a value of the path-angle results which is to high. 
A new fast variable is examined in Ref 43 which picks path- 
angle values in between these two values, and which may be 
used as initial conditions for the problem at hand. 


Chapter XI 

IMPLEMENTATION AND CONCLUSIONS 
11.1 IMPLEMENTATION 

Before the scheme may be used on a real aircraft there are 
some important simplifications and restrictions which have 
been applied in the interest of reducing the initial work- 
load which must be accounted for. 

First, the weight variation of the aircraft must be included 
in the modelling as a substantial percentage of the total 
weight may be used up during a mission. This is perhaps the 
easiest or at least the most straight-forward problem; the 
required action is to increase the order of the system, i.e. 
the mass is added as another variable and the resulting 
boundary conditions are simply that the initial mass is 
known, initial mass multiplier is unknown, and the final 
mass is unknown resulting in the mass multiplier being zero 
at the final time. 

Fuel optimization is a problem which will no doubt be of 
interest, with different combinations of fuel and range 
being optmized. Problems can occur here with a non-convex 
hodograph, i.e. leading to the possibility of chattering 
controls, in this case the throttle. 
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Other problems of the real world which have not been 
addressed are variations in atmospheric conditions i.e. 
winds aloft and non-standard temperature distribution 
against altitude. Possibly these could be dealt with by 
anlaysing the effect of small perturbations, finding an 
approximation to the first order changes in the variables 
which are stored on-board and using simple linear 
corrections. Certainly this is the simplest way of tackling 
such difficulties and it would be interesting to examine how 
effective this approach would be. 

Another problem of interest is that of variable 
configuration, i.e. the effect on the guidance scheme of 
changes in the aircrafts characteristics due to battle 
damage, releasing external stores, etc. 

The biggest problem that must be looked at is the extension 
to 3-D, discussed in the last chapter. 


11.2 CONCLUSIONS 

The' numerical results bear out the following conclusions; 
first, that all trajectories which fair into the high-speed 
point consist of a rapid transition onto a reference or 
skeletal path if they do not originate on it- Secondly, the 
linear-feedback scheme proposed is able to control the 
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aircraft so that it closely follows the appropriate neighbor 
of the nominal path for large perturbations of initial 
conditions. 

11.3 FUTURE WORK 

A 3-D eKtension of the computational scheme is of interest 
in which there are two dominant states, i.e. heading-to-go 
in addition to energy-to-go. As a result, families of 
optimal paths which fair into the dash-point will be needed, 
and the feedback coefficients will be functions of two 
variables (represented via a spline lattice) instead of one. 
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Appwi<tix A 


~)p pOOi^ 


•..i 

Qyi\Uf^ 


Cdo - 0.02&Z + krcun(S0(M-1.0})(l.CHO.3Sexp(>4.S(K-l.S)^H0.0IZ/:r) 
+ 0.08«p(-55(M-l.l)^l + 0.0096«sp(-Z0CM-l.35)^) 

+ 0.003«*p(-20(M-I.6)^) 

Cdel - (0.34O.2026srccan(50CK-1.23)>arctan{S0(Z.25-M»(0.39M-0.475) 

+ 0.07S + O.GSexp(-130(K-O.98S)^) + 0.4(0.S+arcc«n(50(M-2.25)) 

Claw “ 0.82 + (0.72/7r)wctaa(50(0.9-M>)+ 

(1 . 23^0 . 6MH0 . 540 . 2026arctan(SO (H.0 .9})arecan(50(2.0S-H» 

Thruac(M,h) " 

(0.5+<l/jr )«rci:ao<40(M-XM2)HH2-Hl) 4- HI 4- 

(0.5+<2/ jr^)«rcCin(40(M-aU))«rctaxi(60Cai2-«))) (H2-HI/SM2-JH1) (H-XMl) 


2M1,1Q{2,H1,112 are funccions of altitude: 

XMl - (3.84(exp(0.X65((h4-l.74))>) - 4.82) 
nC - 0.0156h^ 4- 2.83K 4-1.1 
HI - (fl.gl 4- f2.B2)£3 (41000) 

H2 “ (fU.gl44E22.g2)(0.54-(l/7r)arctan(40(0.9l-h»)40405 

fl - -2.43h^ - l.S9h 4- 0.974 

f2 - 2.38h^ - 3.24h + 1.24 

gl - (O.S 4- (l/7r)arccaa(40(0.3-h)) 

g2 ■ (0.5 4- (l/n’)arctaa(40(h-0.3)) 

£11 - 1.35li^ - 1.53h 4- 1.56 
£22 - 3.25h^ - 6.25h + 2.98 
£3 • (0.5 4- (l/7r)arctan(40(0.7S - M)) 

h “ altitude/ 10^ 
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Mach Limit Fairing 

The thrust is multiplied by the factor given by; 
f - (0,5 + (I/tt JarctanC150(2.4 - M) 

Dynamic Pressure Limit Fairing 

The thrust is multiplied by the factor given by: 

f - (0.5 + a/7r)arctan(150(M* - M) 

M* « \/(4000/rtio) / as 


riio " density 

S8 « speed of sound 
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Table 1 i 

Odo Da-ta 


1 Mach Number 

1 Cdo 

1 

] 0.00 

1 0.01950 

1 

1 0.50 

1 0.01950 

1 

1 0.80 

1 0.01950 

1 

1 0.88 

1 0.02097 

1 

1 0.90 

i 0.02134 

1 

1 1.00 

1 0.03533 

1 

1 1.10 

1 0 . 04095 

1 

1 1.20 

1 0.04656 

1 

1 1.30 

1 0.04570 

1 

( 1.40 

1 0.04950 

1 

1 1.50 

1 0.04934 

1 

1 1.60 

1 0 . 04918 

1 

1 1.70 

1 0 . 04744 

1 

1 1.80 

1 0.04570 

1 

I 1.90 

1 0.04450 

1 

I 2.00 

I 0.04330 

1 

1 2.10 

1 0.04166 

I 

1 2.20 

1 0.04001 

1 

1 2.30 

1 0.03801 

] 

1 2.50 

1 0.03451 

1 
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Table 2 i 

0dcl2 Data 


1 Mach Number 

i CdC12 

1 

1 0.00 

1 0.07500 

1 

1 0.40 

1 0.07500 

I 

1 0.60 * 

1 0.07500 

1 

1 0.77 

1 0.07500 

i 

I 0.80 

1 0.07500 

1 

1 0.90 

1 0.10000 

1 

I 1.00 

1 0.12500 

1 

1 1.10 

1 0.07500 

1 

1 1.20 

) 0.10000 

1 

1 1.40 

1 0.15000 

1 

1 1.60 

1 0.22500 

1 

r 1.80 

1 0.30000 

1 

1 2.00 

1 0.38750 

1 

1 2.15 

1 0.45000 

1 

( 2.20 

1 0.47500 

1 

1 2.25 

1 0.47500 

1 

1 2.40 

I 0.47500 

1 



Table 3 


Clmax Data 


,<r> 


1 

Hach Number | 

1 Clmax 

i 

1 

0.00 

1 1.180 

1 

1 

0.40 

1.180 

1 

1 

0.60 

1.180 


1 

0.80 

1.160 

1 

1 

1.00 

1.080 

1 

1 

1.20 

0.930 

( 

1 

1.40 

0.810 

1 

1 

1.60 

0.700 

1 

i 

1.80 

0.630 

1 

! 

2.00 

0.570 

1 

( 

2.20 

0.500 

1 

! 

2.40 

0.460 

1 

1 

2.50 

0.460 

! 



Hach Number 


\<i£i 

) |0 

0.2 

0.4 

0.6 

O.B 

1.0 

1-2 

1.4 

1.6 1.8 

2.0 

2.2 

2.4 

2.6 

3.0] 

A 

j 

3B5BO 

41110 

44270 

>19630 

54390 

62170 

65600 

65600 

65600 65600 

65600 

65600 

59040 

29520 

01 

1 


|10 

1 

33070 

33070 

33070 

37160 

42370 

49950 

54230 

59930 

59930 59930 

59930 

59930 

53940 

29520 

1 

0| 

1 

0| 

1 

0| 

1 


1 

|20 

1 

23240 

23240 

23240 

26830 

30900 

37630 

43520 

48070 

54770 54770 

54770 

54770 

49290 

24645 

3* 

2 

1 

[30 

1 

140 

1 

150 

1 

160 

1 

|70 

1 

15540 

1S540 

15540 

17950 

21210 

26450 

32120 

38070 

42600 48010 

52430 

55810 

50230 

25115 

Irt 

ft 

12000 

12000 

12000 

12000 

14400 

18400 

22400 

26800 

30200 34000 

37200 

40000 

36000 

18000 

1 

0| 

1 

0! 

1 

01 

1 

? 

ET 

6648 

6648 

6648 

6648 

7978 

10230 

12800 

15110 

17220 20150 

21960 

24020 

225B0 

11290 
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5660 

5660 

5660 

5660 

5660 

6540 

8200 

9800 11200 

12760 
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13320 
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cr 

lA 

4400 
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4400 5400 

6400 
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6930 
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1 
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1 

01 

1 

01 
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2200 

2200 
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3600 
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Fig 39 


Path- Angle Gain vs Log VS) 
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Fig 48 Interceptor Path without y as Running Variable 




